function p = params

p.Tmax = 50;
p.Tmax = 1e4;

p.N = 300;
p.dz = 1 / (p.N - 1);

p.delta = 0.01;
p.beta = 0.85;
p.epsilon_0 = 0.00;
p.Ma_prime = (1 + p.epsilon_0) * 105;
p.Ma = p.Ma_prime / p.delta / p.beta / (1 - p.beta);